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PACS. 61.20.Ja - Computer simulation of liquid structure. 
PACS. 61.20.Ne - Structure of simple liquids. 
PACS. 64.70.Pf - Glass transitions. 



^ ' Abstract. - Molecular dynamics simulations of a Lennard- Jones binary mixture confined in 

e, a disordered array of soft spheres are presented. The single particle dynamical behavior of the 
glass former is examined upon supercooling. Predictions of mode coupling theory are satisfied 

i-pj I by the confined liquid. Estimates of the crossover temperature are obtained by power law fit 

r^ . to the diffusion coefficients and relaxation times of the late a region. The b exponent of the 

Q ' von Schweidler law is also evaluated. Similarly to the bulk, diff'erent values of the exponent 7 

O , are extracted from the power law fit to the diffusion coefficients and relaxation times. 

(N 
> 



The glass transition scenario in confined liquids is a field of rapidly growing interest because 
vQ ' of its connection with relevant problems in technology and biology H . One of the momentous 
f — . . questions in this field is whether and to what extent theories successfully used for investigating 
^^ ' the glass transition in the bulk retain their validity upon confinement. Experiments M-M show 
a rather diversified phenomenology. Although in most cases the process of vitrification of a 
glass former upon supercooling is not suppressed by confinement, substantial modifications 
with respect to the bulk behavior are almost always observed 

In the theoretical framework of the glass transition of bulk liquids a relevant role is played 

by Mode Couphng Theory (MCT) [|lO|. This theory is able to predict the approach of the 

''O ' supercooled liquid [ pd| that is undergoing a calorimetric glass transition to a temperature Tc 

G . that signals the ideal crossover from a regime where the breaking and the formation of cages 

is the mechanism that ensures the system to be ergotic to a regime where cages are frozen 

and only hopping processes restore ergodicity. 

The approach to the glass transition of binary mixtures in the bulk phase has been the 
^^ . subject of numerous studies, see for example 1 12 - 14| . In particular a numerical test of the MCT 



^H ' predictions has been successfully carried out for a Lennard Jones binary mixture (LJBM) Ml 



Real porous solids are often disordered, they are made of an interconnected network of voids 
of various size and shapes. In order to investigate the influence of a disordered microstructure 
on the glass transition scenario we assume for the solid the same model used in the past both 
in computer simulation and in statistical mechanical theories to study the modifications of 
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the phase diagram in confined fluids: the sohd is treated as a collection of particles forming a 
rigid disordered structure pq-Eq]. In this way the fluid is embedded in an off- lattice matrix 
of rigid "obstacles" . This model is in particular representative of amorphous materials like 
silica xerogels. 

In this letter we demonstrate that the main predictions of MCT are satisfied for a glass 
former liquid confined in this disordered medium. 

We introduced into a rigid disordered matrix of 16 soft spheres the LJBM of 800 A and 
200 B particles as defined in ref^JnSl. The use of a binary mixture is strictly necessary in 
the bulk to avoid cristallization [[19[. Confinment is expected to inhibit cristallization also in 
a one component fluid instead. We did however choose to use the LJBM of ref. [Ea for a 
comparison with the bulk. The parameters of the LJBM are eAA = 1, o'aa ~ 1, ^bb — 0.5, 
UBB = 0.88, eAB = 1-5, and aAB = 0.8. In the following energy will be given in units of eaa, 
temperature in units of caa/Kb, length in units <jaa and time in units of (mcr^^/(48e^^))^'^. 
The Lennard-Jones particles interact with the soft spheres with a repulsive potential: 



Vir)=Ae[^) (1) 

with esA — 0.32, asA — 3, €sb = 0.22, asB — 2.94. Parameters have been chosen in order to 
obtain a strongly confined system where most of the L J particles are interfacial. The potential 
is truncated at r|- — 2.5aij. We conducted the simulations in the microcanonical ensemble in 
a cubic box of length 12.6 with periodic boundary conditions. The equations of motion were 
solved by the velocity Verlet algorithm. 

The system was equilibrated at different reduced temperatures via a velocity rescaling 
procedure starting from T ^ 5.0. It was then progressively cooled and equilibrated at the 
following temperatures T = 4.0, 3.0, 2.0, 1.0. Below T = 1.0, after equilibration, production 
runs were done for the following temperatures: T = 0.80, 0.58, 0.48, 0.43, 0.41, 0.39 and 0.37. 
The timestep used for T > 1 was 0.01 and for T < 1 was 0.02. For the lowest temperature 
investigated a production run of ten million timesteps was performed. We verified that the 
results presented in this letter do not depend on the specific choice of the disordered matrix 
by running MD simulations for other two different configurations of the disordered matrix of 
soft spheres. Both thermodynamics and dynamics appeared the same for the three systems, 
in particular deviations among the relaxation times extracted at a given temperature from 
the density correlators of different configurations are within 2%. 

A snapshot of the system under investigation is shown for a low temperature in Fig. 1. 
From the picture it is evident that we are dealing with a strongly confined system. 

No phase separation of A and B particles is detected during the cooling process through 
the calculation of the pair correlation functions (not shown) Ig^l ■ 

In the following we concentrate on the single particle dynamics of A and B particles. In 
Fig. 2 the Mean Square Displacement (MSD) for the A particles is shown for all the temper- 
atures investigated. It displays the typical behavior of a system approaching the crossover 
temperature of MCT. Upon lowering the temperature it appears an intermediate time region 
where the MSD flattens due to the well known cage effect. Similar behavior is found for the 
B particles (not shown). From the slope of the MSD, once the particle enters the Brownian 
diffusive regime, it is possible to extract the diffusion coefficient which is predicted by MCT 
to have a power law behavior: 

D^{T~ TcV (2) 

In the inset of Fig. 2 we show the diffusion coefficients extracted at the different temperatures 
from the MSD both for A and B particles together with the fit to Eq. 0. For A particles 
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Tc — 0.343 and 7 = 1.60 and for B particles Tc = 0.343 and 7 — 1.69 are found. We obtain 
the same value of Tc for both species, as predicted by MCT, but slightly different values of 
7. This discrepancy was also found in the bulk pa. 

In Fig. 3 we show the density-density self correlation function Fs{Q,t) at the peak of the 
AA structure factor, which for this system is Qctaa — 7.07 as shown in the inset of Fig. 4. We 
calculated Fs{Q,t) also for B particles at the peak of the BB structure factor Qctbb ~ 5.90 
(not shown). As already displayed by the MSD for higher temperatures the function does show 
a single relaxation timescale while upon cooling we observe the shouldering of the relaxation 
law that corresponds to the onset of the slow a-relaxation regime in the system. The late part 
of the a relaxation region is found in most glass formers to have the analytical shape of the 
Kohlrausch- William- Watts, KWW, function: 

Mt) = fQe-'^'^^^' (3) 

where /g is the height of the plateau and r is the late a relaxation time. The fit to this 
formula are also shown in Fig. 3. The agreement is very satisfactory. In the lower inset we 
report the f3 and fg values extracted for the A and B particles from the fit. In the upper 
inset the r values are shown together with the fit to the MCT power law: 

T^{T-T,)-^ (4) 

For A particles Tc = 0.343 and 7 = 2.90 and for B particles Tc = 0.343 and 7 = 2.89 are 
extracted. As found in the bulk the values of Tc predicted from D and r are similar while 
the 7 are different. 

We now move to the so called von Schweidler test. MCT predicts that the region of 
departure from the plateau of the density correlator in the Q,t space, the late /3 relaxation 
region, behaves according to the following power law: 

Mt) - fQ - hQit/r)" (5) 

where r is a characteristic time of the system, b a system dependent exponent, /q is the height 
of the plateau at the temperature Tc and hq an amplitude that does not depend on time. In 
Fig. 4 we show the Fs{Q, t) of Fig. 3 now rescaled according to the r values extracted from 
the fit to Eq. y. We note that all the curves collapse into a single master curve as predicted 
by MCT. The deviations from the master curve appear to be mostly due to overshots clearly 
visible in our correlators and also present in the bulk that tend to mask the MCT predicted 
behavior. This kind of oscillations appear to be more marked upon confinement. In the same 
figure the fit of the master curve to Eq. g is also reported. We obtain b — 0.355 and /^ = 0.72 
for A particles, b — 0.35 and f^ = 0.785 for B particles. 

MCT also predicts that the region of approach to the plateau is given by another power 
law characterized by an exponent a, also called the critical exponent, a cannot be directly 
extracted by our data due to the oscillations that start appearing in the correlator upon 
cooling, as also evident in Fig. 3. 

With the values of b and the 7 extracted by Eq. H and by means of the MCT relationship: 

^-i + Yb (^) 

we obtain a = 0.335 for the A particles and a = 0.342 for the B particles. The values of 7, 
a and b are well in the range predicted by the theory. The exponents a and b both allow an 
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evaluation of the parameter A of MCT, also known as the exponent parameter, 1/2 < A < 1, 
by means of the MCT formula: 

^ [T{l-a)f ^ [r(l + b)f 

r(l-2a) r(i + 2&) ^ ' 

where r(a;) is the T function. 

The system dependent parameter A governs the behaviour of /g close to /q below Tc- 
We have a difference of 22% for the values of A extracted from b and a for both A and B 
species which is larger with respect to the value 5% found for the bulk [^. In our case we 
observe that, due to the inevitable oscillations of the correlators, the fit of the von Schweidler 
law that allow the determination of the exponent b is affected by a larger uncertainty with 
respect to the bulk case. Therefore the best estimate we can obtain for the parameter A of 
this system is: A = 0.8 ±0.1 where we are attributing to A the maximum error coming form 
the semidispersion of the values obtained. We note that for a more precise determination of 
the exponent b also the term next to the leading order has to be taken into account in eq.^ 
for the fit and at the same time several Q values should be tested simultaneously |2^. This 
analysis is in progress and will be reported in a subsequent paper [pO| . 

In conclusion we performed a MD simulation on a LJBM confined in a disordered matrix 
of soft spheres. Our results provide insight into the microscopic behavior of this confined fluid 
when cooled. Since no cristalline equilibrium solid can be produced in the confined matrix we 
are not able to individuate the freezing temperature below which the liquid can be considered 



supercooled as it is done in the bulk [14| . None-the-less most of liquids start displaying the 
shouldering of the relaxation laws below the freezing temperature, we can therefore reasonably 
consider the fluid supercooled for temperatures below T = 0.5. We have shown that MCT is 
able to rationalize the dynamical behavior of the system in spite of the fact that this theory 
is formulated for bulk liquids. In particular both the diffusion coefficients extracted from the 
slope of the MSD and the relaxation times extracted from the fit to the KWW function to the 
late part of the a relaxation region of the intermediate scattering function allow to estimate the 
crossover temperature Tc = 0.343 of the confined LJBM. The rescaled intermediate scattering 
functions define a master curve that can be fitted to the von Schweidler law from which the b 
exponent is extracted. From the values of a and b the exponent parameter A is also evaluated. 

An MCT test with MD on confined liquids has been so far carried out to the best of our 
knowledge only for water confined in a silica nanopore [2^ . In that case due to the different 
kind of confining medium and to the different liquid the MCT behaviour could be extracted 
only with a layer analysis of the density correlators. For the system presented here this is 
not necessary. Besides the obvious differences in the confining geometry and liquid another 
important difference between the two systems is that the LJBM is in a situation of strong 
confinement. In fact given the paramenter of the soft spheres only roughly 1/3 of the total 
volume is accessible for the mixture and only few layers of LJBM can reside among the spheres 
as it is also evident from Fig. 1. In the case of water many layers can reside in a nanopore. 

Although our system is not directly comparable to any experimental system, simple fiuids 
adsorbed in silica xerogels are likely to display similar behavior p6|-p8[ . Experimental MCT 
tests on such systems would therefore be very valuable to assess the validity of our findings. 
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Fig. 1 - Snapshot of the system simulated at T = 0.4. Larger spheres represent the confining system, 
dark-gray particles the A type, light-gray the B type particles of the LJBM. 
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Fig. 2 - Mean square displacement of A particles for T = 0.80, 0.58, 0.48, 0.43, 0.41, 0.39, 0.37. 
Curves on the top correspond to higher temperatures. In the inset we show the diffusion coefficient D 
as a function of temperature for A (squares) and B (circles) particles together with the fit (continuous 
lines) to the power law of Eq. eI For A particles Tc = 0.343 and 7 = 1.60, for B particles Tc = 0.343 
and 7 — 1.69. 



EUROPHYSICS LETTERS 




10 

units) 



Fig. 3 - Self part of the intermediate scattering function at the peak of the structure factor Qmax = 
7.07a^]^ for A particles (continuous line) for the same temperatures as in Fig. ti. Curves on the top 
correspond to lower temperatures. The dashed curves are the KWW fits. In the lower inset we show 
the non-ergodicity parameter /q (filled symbols) and the stretching exponent /3 (empty symbols) for 
A (squares) and B (circles) particles, as extracted from the fits to the KWW. In the upper inset 
the power law fit (continuous lines) for the relaxation time is shown for A (squares) and B (circles) 
particles. For A particles Tc = 0.343 and 7 = 2.90, for B particles Tc = 0.343 and 7 = 2.89. 
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Fig. 4 - In the main frame continuous lines are the time-rescaled intermediate scattering functions 
of Fig. M. The dashed Une is the fit to the master curve according to the von Schweidler law. We 
obtain b = 0.355 and fg — 0.72 for A particles, b = 0.35 and /q = 0.785 for B particles (not shown). 
In the inset the static structure factor of A particles is shown for all the temperatures investigated. 
For the sake of clarity lower temperatures have been upward shifted. The first peak is located at 
Qmax = 7.07cr] 
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